#include "sampling.hh"
#include <cstdlib>
#include <cmath>

using namespace cgmath;

vector_3d lambert_sample ()
{
  double x, y, z;
  do
    {
      x = 2.0 * (rand() / (double)RAND_MAX) - 1.0;
      y = 2.0 * (rand() / (double)RAND_MAX) - 1.0;
    }
  while (x*x + y*y > 1.0);
  z = sqrt (1.0 - x*x - y*y);
  return vec (x, y, z);
}

vector_3d phong_sample (double ns)
{
  double r1 = rand() / (double) RAND_MAX;
  double r2 = rand() / (double) RAND_MAX;

  double x = cos(2*M_PI*r1)*sqrt(1-pow(r2, (2.0/(ns+1))));
  double y = sin(2*M_PI*r2)*sqrt(1-pow(r2, (2.0/(ns+1))));
  double z = pow(r2, 1.0/(ns+1));

  return vec (x, y, z);
}

// We pick points at the square [-1, 1] x [-1, 1]
// uniformly distributed until we find one that is inside the unit
// disk.
vector_3d lens_sample()
{
  vector_3d p;
  double sx, sy;
  while (true)
    {
      sx = 2 * (double)rand() / (double)RAND_MAX - 1.0;
      sy = 2 * (double)rand() / (double)RAND_MAX - 1.0;

      if (sx * sx + sy * sy <= 1)
	{
	  break;
	}
    }

  p(0) = sx;
  p(1) = sy;
  p(2) = 0.0;

  return p;
}

vector_3d mirror_sample (vector_3d normal, vector_3d out_dir)
{
  matrix_3d R = rotate (M_PI, normal);
  return R * out_dir;
}
